#####################################################################
# Outcomes
#####################################################################

rm(list=ls())

sink("~/Dropbox/Geographic Natural Experiments/02_analyses/05_match_balance.txt")

library(Hmisc)
library(exactRankTests)
library(ggplot2)
library(reshape2)
library(xtable)
library(stargazer)
library(sensitivitymw)
library(sensitivitymv)

##############################################
# Read data and prepare basic output
##############################################

load("~/Dropbox/Geographic Natural Experiments/02_analyses/01_match_1000m.RData")
load("~/Dropbox/Geographic Natural Experiments/02_analyses/02_match_city_1000m.RData")
load("~/Dropbox/Geographic Natural Experiments/02_analyses/03_match_state_1000m.RData")
load("~/Dropbox/Geographic Natural Experiments/02_analyses/04_match_country_1000m.RData")
names(d)

#########################
# Explore data
#########################

# observations
table(d$t_ind)
table(d_match$t_ind)
table(d_match_state$t_ind)
table(d_match_city$t_ind)
table(d_match_country$t_ind)

# median vs mean
mean(d$age)
median(d$age)

#############################################
# Standardized differences before matching
#############################################

num1 = abs(mean(d$hprice_mean[d$t_ind==1])- mean(d$hprice_mean[d$t_ind==0]))
den1 = sqrt(  ( (sd(d$hprice_mean[d$t_ind==1]))^2 + (sd(d$hprice_mean[d$t_ind==0]))^2 )/2  )
sd1 = num1/den1
sd1_r = round(sd1,2)
sd1_r

num2 = abs(mean(d$turnout04[d$t_ind==1])- mean(d$turnout04[d$t_ind==0]))
den2 = sqrt(  ( (sd(d$turnout04[d$t_ind==1]))^2 + (sd(d$turnout04[d$t_ind==0]))^2 )/2  )
sd2 = num2/den2
sd2_r = round(sd2,2)
sd2_r

num3 = abs(mean(d$turnout06[d$t_ind==1])- mean(d$turnout06[d$t_ind==0]))
den3 = sqrt(  ( (sd(d$turnout06[d$t_ind==1]))^2 + (sd(d$turnout06[d$t_ind==0]))^2 )/2  )
sd3 = num3/den3
sd3_r = round(sd3,2)
sd3_r

num4 = abs(mean(d$male[d$t_ind==1])- mean(d$male[d$t_ind==0]))
den4 = sqrt(  ( (sd(d$male[d$t_ind==1]))^2 + (sd(d$male[d$t_ind==0]))^2 )/2  )
sd4 = num4/den4
sd4_r = round(sd4,2)
sd4_r

num5 = abs(mean(d$age[d$t_ind==1])- mean(d$age[d$t_ind==0]))
den5 = sqrt(  ( (sd(d$age[d$t_ind==1]))^2 + (sd(d$age[d$t_ind==0]))^2 )/2  )
sd5 = num5/den5
sd5_r = round(sd5,2)
sd5_r

##################################################
# Standardized differences after regular matching
##################################################

num1b = abs(mean(d_match$hprice_mean[d_match$t_ind==1])- mean(d_match$hprice_mean[d_match$t_ind==0]))
sd1b = num1b/den1
sd1b_r = round(sd1b,2)
sd1b_r

num2b = abs(mean(d_match$turnout04[d_match$t_ind==1])- mean(d_match$turnout04[d_match$t_ind==0]))
sd2b = num2b/den2
sd2b_r = round(sd2b,2)
sd2b_r

num3b = abs(mean(d_match$turnout06[d_match$t_ind==1])- mean(d_match$turnout06[d_match$t_ind==0]))
sd3b = num3b/den3
sd3b_r = round(sd3b,2)
sd3b_r

num4b = abs(mean(d_match$male[d_match$t_ind==1])- mean(d_match$male[d_match$t_ind==0]))
sd4b = num4b/den4
sd4b_r = round(sd4b,2)
sd4b_r

num5b = abs(mean(d_match$age[d_match$t_ind==1])- mean(d_match$age[d_match$t_ind==0]))
sd5b = num5b/den5
sd5b_r = round(sd5b,2)
sd5b_r

#########################################################
# Standardized differences after template matching: city
#########################################################

num1c = abs(mean(d_match_city$hprice_mean[d_match_city$t_ind==1])- mean(d_match_city$hprice_mean[d_match_city$t_ind==0]))
sd1c = num1c/den1
sd1c_r = round(sd1c,2)
sd1c_r

num2c = abs(mean(d_match_city$turnout04[d_match_city$t_ind==1])- mean(d_match_city$turnout04[d_match_city$t_ind==0]))
sd2c = num2c/den2
sd2c_r = round(sd2c,2)
sd2c_r

num3c = abs(mean(d_match_city$turnout06[d_match_city$t_ind==1])- mean(d_match_city$turnout06[d_match_city$t_ind==0]))
sd3c = num3c/den3
sd3c_r = round(sd3c,2)
sd3c_r

num4c = abs(mean(d_match_city$male[d_match_city$t_ind==1])- mean(d_match_city$male[d_match_city$t_ind==0]))
sd4c = num4c/den4
sd4c_r = round(sd4c,2)
sd4c_r

num5c = abs(mean(d_match_city$age[d_match_city$t_ind==1])- mean(d_match_city$age[d_match_city$t_ind==0]))
sd5c = num5c/den5
sd5c_r = round(sd5c,2)
sd5c_r

###########################################################
# Standardized differences after template matching: state
###########################################################

num1d = abs(mean(d_match_state$hprice_mean[d_match_state$t_ind==1])- mean(d_match_state$hprice_mean[d_match_state$t_ind==0]))
sd1d = num1d/den1
sd1d_r = round(sd1d,2)
sd1d_r

num2d = abs(mean(d_match_state$turnout04[d_match_state$t_ind==1])- mean(d_match_state$turnout04[d_match_state$t_ind==0]))
sd2d = num2d/den2
sd2d_r = round(sd2d,2)
sd2d_r

num3d = abs(mean(d_match_state$turnout06[d_match_state$t_ind==1])- mean(d_match_state$turnout06[d_match_state$t_ind==0]))
sd3d = num3d/den3
sd3d_r = round(sd3d,2)
sd3d_r

num4d = abs(mean(d_match_state$male[d_match_state$t_ind==1])- mean(d_match_state$male[d_match_state$t_ind==0]))
sd4d = num4d/den4
sd4d_r = round(sd4d,2)
sd4d_r

num5d = abs(mean(d_match_state$age[d_match_state$t_ind==1])- mean(d_match_state$age[d_match_state$t_ind==0]))
sd5d = num5d/den5
sd5d_r = round(sd5d,2)
sd5d_r

###########################################################
# Standardized differences after template matching: country
###########################################################

num1e = abs(mean(d_match_country$hprice_mean[d_match_country$t_ind==1])- mean(d_match_country$hprice_mean[d_match_country$t_ind==0]))
sd1e = num1e/den1
sd1e_r = round(sd1e,2)
sd1e_r

num2e = abs(mean(d_match_country$turnout04[d_match_country$t_ind==1])- mean(d_match_country$turnout04[d_match_country$t_ind==0]))
sd2e = num2e/den2
sd2e_r = round(sd2e,2)
sd2e_r

num3e = abs(mean(d_match_country$turnout06[d_match_country$t_ind==1])- mean(d_match_country$turnout06[d_match_country$t_ind==0]))
sd3e = num3e/den3
sd3e_r = round(sd3e,2)
sd3e_r

num4e = abs(mean(d_match_country$male[d_match_country$t_ind==1])- mean(d_match_country$male[d_match_country$t_ind==0]))
sd4e = num4e/den4
sd4e_r = round(sd4e,2)
sd4e_r

num5e = abs(mean(d_match_country$age[d_match_country$t_ind==1])- mean(d_match_country$age[d_match_country$t_ind==0]))
sd5e = num5e/den5
sd5e_r = round(sd5e,2)
sd5e_r

##############################################
# Balance 
##############################################

# Before matching 
treated_before_house = mean(d$hprice_mean[d$t_ind==1])
control_before_house = mean(d$hprice_mean[d$t_ind==0])
treated_before_turnout04 = mean(d$turnout04[d$t_ind==1])
control_before_turnout04 = mean(d$turnout04[d$t_ind==0])
treated_before_turnout06 = mean(d$turnout06[d$t_ind==1])
control_before_turnout06 = mean(d$turnout06[d$t_ind==0])
treated_before_male = mean(d$male[d$t_ind==1])
control_before_male = mean(d$male[d$t_ind==0])
treated_before_age = mean(d$age[d$t_ind==1])
control_before_age = mean(d$age[d$t_ind==0])

# After matching 
treated_after_house = mean(d_match$hprice_mean[d_match$t_ind==1])
control_after_house = mean(d_match$hprice_mean[d_match$t_ind==0])
treated_after_turnout04 = mean(d_match$turnout04[d_match$t_ind==1])
control_after_turnout04 = mean(d_match$turnout04[d_match$t_ind==0])
treated_after_turnout06 = mean(d_match$turnout06[d_match$t_ind==1])
control_after_turnout06 = mean(d_match$turnout06[d_match$t_ind==0])
treated_after_male = mean(d_match$male[d_match$t_ind==1])
control_after_male = mean(d_match$male[d_match$t_ind==0])
treated_after_age = mean(d_match$age[d_match$t_ind==1])
control_after_age = mean(d_match$age[d_match$t_ind==0])

# After matching city
treated_after_city_house = mean(d_match_city$hprice_mean[d_match_city$t_ind==1])
control_after_city_house = mean(d_match_city$hprice_mean[d_match_city$t_ind==0])
treated_after_city_turnout04 = mean(d_match_city$turnout04[d_match_city$t_ind==1])
control_after_city_turnout04 = mean(d_match_city$turnout04[d_match_city$t_ind==0])
treated_after_city_turnout06 = mean(d_match_city$turnout06[d_match_city$t_ind==1])
control_after_city_turnout06 = mean(d_match_city$turnout06[d_match_city$t_ind==0])
treated_after_city_male = mean(d_match_city$male[d_match_city$t_ind==1])
control_after_city_male = mean(d_match_city$male[d_match_city$t_ind==0])
treated_after_city_age = mean(d_match_city$age[d_match_city$t_ind==1])
control_after_city_age = mean(d_match_city$age[d_match_city$t_ind==0])

# After matching state
treated_after_state_house = mean(d_match_state$hprice_mean[d_match_state$t_ind==1])
control_after_state_house = mean(d_match_state$hprice_mean[d_match_state$t_ind==0])
treated_after_state_turnout04 = mean(d_match_state$turnout04[d_match_state$t_ind==1])
control_after_state_turnout04 = mean(d_match_state$turnout04[d_match_state$t_ind==0])
treated_after_state_turnout06 = mean(d_match_state$turnout06[d_match_state$t_ind==1])
control_after_state_turnout06 = mean(d_match_state$turnout06[d_match_state$t_ind==0])
treated_after_state_male = mean(d_match_state$male[d_match_state$t_ind==1])
control_after_state_male = mean(d_match_state$male[d_match_state$t_ind==0])
treated_after_state_age = mean(d_match_state$age[d_match_state$t_ind==1])
control_after_state_age = mean(d_match_state$age[d_match_state$t_ind==0])

# After matching country
treated_after_country_house = mean(d_match_country$hprice_mean[d_match_country$t_ind==1])
control_after_country_house = mean(d_match_country$hprice_mean[d_match_country$t_ind==0])
treated_after_country_turnout04 = mean(d_match_country$turnout04[d_match_country$t_ind==1])
control_after_country_turnout04 = mean(d_match_country$turnout04[d_match_country$t_ind==0])
treated_after_country_turnout06 = mean(d_match_country$turnout06[d_match_country$t_ind==1])
control_after_country_turnout06 = mean(d_match_country$turnout06[d_match_country$t_ind==0])
treated_after_country_male = mean(d_match_country$male[d_match_country$t_ind==1])
control_after_country_male = mean(d_match_country$male[d_match_country$t_ind==0])
treated_after_country_age = mean(d_match_country$age[d_match_country$t_ind==1])
control_after_country_age = mean(d_match_country$age[d_match_country$t_ind==0])

# Mean vectors
mom_targets_city = round(c(152996,0.6998,0.569,0.478,30),2)
mom_targets_state = round(c(162407,0.729,0.509,0.494,36),2)
mom_targets_country = round(c(214546,0.638,0.478,0.491,35.3),2)
before_treated = round(c(treated_before_house,treated_before_turnout04,treated_before_turnout06,treated_before_male,treated_before_age),2)
before_control = round(c(control_before_house,control_before_turnout04,control_before_turnout06,control_before_male,control_before_age),2)
after_treated = round(c(treated_after_house,treated_after_turnout04,treated_after_turnout06,treated_after_male,treated_after_age),2)
after_control = round(c(control_after_house,control_after_turnout04,control_after_turnout06,control_after_male,control_after_age),2)
after_treated_city = round(c(treated_after_city_house,treated_after_city_turnout04,treated_after_city_turnout06,treated_after_city_male,treated_after_city_age),2)
after_control_city = round(c(control_after_city_house,control_after_city_turnout04,control_after_city_turnout06,control_after_city_male,control_after_city_age),2)
after_treated_state = round(c(treated_after_state_house,treated_after_state_turnout04,treated_after_state_turnout06,treated_after_state_male,treated_after_state_age),2)
after_control_state = round(c(control_after_state_house,control_after_state_turnout04,control_after_state_turnout06,control_after_state_male,control_after_state_age),2)
after_treated_country = round(c(treated_after_country_house,treated_after_country_turnout04,treated_after_country_turnout06,treated_after_country_male,treated_after_country_age),2)
after_control_country = round(c(control_after_country_house,control_after_country_turnout04,control_after_country_turnout06,control_after_country_male,control_after_country_age),2)

# standardized differences
sd_before = round(c(sd1_r,sd2_r,sd3_r,sd4_r,sd5_r),2)
sd_regular = round(c(sd1b_r,sd2b_r,sd3b_r,sd4b_r,sd5b_r),2)
sd_city = round(c(sd1c_r,sd2c_r,sd3c_r,sd4c_r,sd5c_r),2)
sd_state = round(c(sd1d_r,sd2d_r,sd3d_r,sd4d_r,sd5d_r),2)
sd_country = round(c(sd1e_r,sd2e_r,sd3e_r,sd4e_r,sd5e_r),2)

# Table 1
colnames = c("Housing prices","Turnout 2004","Turnout 2006","Male","Age")
table1 = data.frame(colnames,before_treated,before_control)
colnames(table1) <- c("Covariates", "Treated","Control")
table1

# Table 2
table2 = data.frame(colnames,after_treated,after_control,mom_targets_city,mom_targets_state,mom_targets_country)
colnames(table2) <- c("Covariates", "Treated","Control","City","State","Country")
table2

# Table 3
table3 = data.frame(colnames,after_treated_state,after_control_state,mom_targets_state)
colnames(table3) <- c("Covariates", "Treated","Control","State")
table3

# Table 4
table4 = data.frame(colnames,after_treated_city,after_control_city,mom_targets_city)
colnames(table4) <- c("Covariates", "Treated","Control","City")
table4

# Table 5
table5 = data.frame(colnames,after_treated_country,after_control_country,mom_targets_country)
colnames(table5) <- c("Covariates", "Treated","Control","Country")
table5

# Table 6
table6 = data.frame(colnames,sd_before)
colnames(table6) <- c("Covariates", "Standardized differences")
table6

# Table 7
table7 = data.frame(colnames,sd_regular)
colnames(table7) <- c("Covariates", "Standardized differences")
table7

# Table 8
table8 = data.frame(colnames,sd_city)
colnames(table8) <- c("Covariates", "Standardized differences")
table8

# Table 9
table9 = data.frame(colnames,sd_state)
colnames(table9) <- c("Covariates", "Standardized differences")
table9

# Table 10
table10 = data.frame(colnames,sd_country)
colnames(table10) <- c("Covariates", "Standardized differences")
table10

# Save table
stargazer(table1, 
          type = "text",
          colnames = TRUE,
          summary = FALSE,
          title="Before matching", 
          digits=2, 
          rownames=FALSE, 
          float = TRUE, 
          float.env = "table", 
          table.placement = "H",
          out="~/Dropbox/Geographic Natural Experiments/03_manuscript/tables/t1_before_matching.tex")

stargazer(table2, 
          type = "text",
          colnames = TRUE,
          summary = FALSE,
          title="After regular matching", 
          digits=2, 
          rownames=FALSE, 
          float = TRUE, 
          float.env = "table", 
          table.placement = "H",
          out="~/Dropbox/Geographic Natural Experiments/03_manuscript/tables/t2_after_matching.tex")

stargazer(table3, 
          type = "text",
          colnames = TRUE,
          summary = FALSE,
          title="After template matching: state", 
          digits=2, 
          rownames=FALSE, 
          float = TRUE, 
          float.env = "table", 
          table.placement = "H",
          out="~/Dropbox/Geographic Natural Experiments/03_manuscript/tables/t3_after_matching_state.tex")

stargazer(table4, 
          type = "text",
          colnames = TRUE,
          summary = FALSE,
          title="After template matching: city", 
          digits=2, 
          rownames=FALSE, 
          float = TRUE, 
          float.env = "table", 
          table.placement = "H",
          out="~/Dropbox/Geographic Natural Experiments/03_manuscript/tables/t4_after_matching_city.tex")

stargazer(table5, 
          type = "text",
          colnames = TRUE,
          summary = FALSE,
          title="After template matching: country", 
          digits=2, 
          rownames=FALSE, 
          float = TRUE, 
          float.env = "table", 
          table.placement = "H",
          out="~/Dropbox/Geographic Natural Experiments/03_manuscript/tables/t5_after_matching_country.tex")

stargazer(table6, 
          type = "text",
          colnames = TRUE,
          summary = FALSE,
          title="Balance before matching", 
          digits=2, 
          rownames=FALSE, 
          float = TRUE, 
          float.env = "table", 
          table.placement = "H",
          out="~/Dropbox/Geographic Natural Experiments/03_manuscript/tables/t6_balance_before.tex")

stargazer(table7, 
          type = "text",
          colnames = TRUE,
          summary = FALSE,
          title="Balance regular matching", 
          digits=2, 
          rownames=FALSE, 
          float = TRUE, 
          float.env = "table", 
          table.placement = "H",
          out="~/Dropbox/Geographic Natural Experiments/03_manuscript/tables/t7_balance_regular.tex")

stargazer(table8, 
          type = "text",
          colnames = TRUE,
          summary = FALSE,
          title="Balance template matching: city", 
          digits=2, 
          rownames=FALSE, 
          float = TRUE, 
          float.env = "table", 
          table.placement = "H",
          out="~/Dropbox/Geographic Natural Experiments/03_manuscript/tables/t8_balance_city.tex")

stargazer(table9, 
          type = "text",
          colnames = TRUE,
          summary = FALSE,
          title="Balance template matching: state", 
          digits=2, 
          rownames=FALSE, 
          float = TRUE, 
          float.env = "table", 
          table.placement = "H",
          out="~/Dropbox/Geographic Natural Experiments/03_manuscript/tables/t9_balance_state.tex")

stargazer(table10, 
          type = "text",
          colnames = TRUE,
          summary = FALSE,
          title="Balance template matching: country", 
          digits=2, 
          rownames=FALSE, 
          float = TRUE, 
          float.env = "table", 
          table.placement = "H",
          out="~/Dropbox/Geographic Natural Experiments/03_manuscript/tables/t10_balance_country.tex")



